#pragma once
#include "feDefs.h"
#include "feMatrix.h"
#include "feVector.h"

#define SOLVER_MAX_ITER 10000
#define BIG_NUMBER 1e20
#define SOLVER_EPS 1e-10

class feVector;
class feMeshFunction;
class feMesh;
class feBordCond;
class feMatrix;
class feFiveDiagSymMatrix;
class feSolver;

class fe1dParabolicTask 
{
	feMeshFunction* Right;

	feMesh  *Mesh;

	double  FirstTimeStep;
	double  CoeffTime;
	int		LayersTime;
	int		LayersPass;

	int  CurLayer;

	feFiveDiagSymMatrix Stiffness;
	feFiveDiagSymMatrix WeightWithGamma;
	feFiveDiagSymMatrix WeightWithDuDt;
	feFiveDiagSymMatrix M;
	feVector			F;
	feVector			*X0;
	feVector			*X1;
	feVector			*X2;

	feSolver			*Solver;

	feBordCond			**Bcs;
	int					NumBc;

private:
	void RunSolver();
	void ProcessBorderConds();

	void SolveNoLayers();
	void Solve2Layer(double iDt);
	void Solve3Layer(double iDt1, double iDt2);

	void ZeroData();
	void GatherStiffness();
	void GatherWeight();
	void GatherRight();

public:
	fe1dParabolicTask();
	bool Init(feMesh *ioMesh, feBordCond **iConds, int iCondsCount, feVector *iX0, double iFirstTimeStep, int iLayersTime, double iStepCoeffTime);
	bool Init(feMesh *ioMesh, feBordCond **iConds, int iCondsCount, feVector *iX0, feMeshFunction *iRight, int iLayersPass, int iNumL, double iScaleCoeff);
	bool SolveNonStationary(feMeshFunction *oFunc);
	bool SolveStationary(feVector *oSolvation, double *oFdotX=0);
	~fe1dParabolicTask();
};
